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Abstract 

The effect of a nonuniform mean flow on (i) the normal modes, (ii) the inflow/outflow 
nonreflecting boundary conditions, and (iii) the sound power are studied. The normal 
modes in an annular duct are computed using a spectral method in combination with 
a shooting method. The swirl causes force imbalance which couples the acoustic and 
vortical modes. The acoustic modes are distinguished from the vortical modes by their 
large pressure and small vorticity content. The mean swirl also produces a Doppler 
shift in frequency. This results in more counter-spinning modes cut-on at a given fre- 
quency than modes spinning with the swirl. Nonreflecting boundary conditions are 
formulated using the normal mode solutions. The inflow/outflow boundary conditions 
are implemented in a linearized Euler scheme and validated by computing the prop- 
agation of acoustic and vortical waves in a duct for a variety of swirling mean flows. 
Numerical results show that the evolution of the vortical disturbances is sensitive to 
the inflow conditions and the details of the wake excitations. All three components of 
the wake velocity must be considered to correctly compute the wake evolution and the 
blade upwash. For high frequencies, the acoustic-vortical mode coupling is weak and a 
conservation equation for the acoustic energy can be derived. Sound power calculations 
show significant mean flow swirl effects, but mode interference effects are small. 
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1 Introduction 


Wakes shed from a fan interact with downstream lying guide vanes and struts to produce 
fan interaction noise. For tonal interaction noise, the harmonic content of the wakes are 
steady in a reference frame moving with the fan but in the absolute reference frame of the 
exit guide vanes, the wakes appear as unsteady disturbances rotating with the blade passing 
frequency of the fan. Our ability to directly compute this unsteady interaction is hindered 
by mathematical and numerical difficulties as well as by the computational intensiveness of 
the problem. A common approach is to limit the size of the computational domain and to 
confine the calculations to a single blade row. The mathematical formulation of the problem 
requires that the inflow and outflow conditions be defined accurately so as to represent both 
the wake harmonics and the upstream/downstream propagating acoustic waves. Moreover, 
causality imposes that a computational domain enclosing the noise source cannot have en- 
ergy propagating into the domain from the inflow/outflow boundaries. This is accomplished 
by imposing nonreflecting inflow/outflow conditions. Since the radiated acoustic energy is 
a very small fraction of the total flow energy, a highly accurate procedure for simulating 
the noise generation mechanism and the propagation of sound must be developed. 

In fan tonal interaction noise, the stator is the isolated blade row. In Figure 1, results 
from a Reynolds Averaged Navier-Stokes calculation for a representative fan design are 
given for the circumferentially averaged mean flow at the inflow plane lying between the 
fan and the stator. In the top figures, the radial profile of the total enthalpy and axial 
velocity are plotted from hub to tip and in the bottom figures the swirl velocity and total 
velocity distribution are plotted from hub to tip. The results show that the mean flow 
variables are nonuniform exhibiting three-dimensional effects with significant swirl velocity. 
The radial variations result mainly from nonuniform blade loading on the fan and viscous 
effects near the hub and tip. An accurate aerodynamic/aeroacoustic computational scheme 
must incorporate upstream and downstream mean flow nonuniformity effects. In particular, 
inflow/outflow conditions must be derived for flows with mean swirl and radial variations. 
This is increasingly important as fan and exit guide vane designs incorporate sweep and lean. 
The main objective of this work is to develop exact nonreflecting inflow/outflow conditions 
for nonuniform mean flows and to implement these conditions for aerodynamic/aeroacoustic 
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calculations in an annular duct geometry. 

Traditionally [1, 2, 3, 4], unsteady linearized Euler algorithms assume uniform mean 
flow conditions at the inflow and outflow boundaries. In this case, it is possible to express 
impinging disturbances as a superposition of acoustic, entropic and vortical modes obeying 
distinct equations [5]. Such approaches have been successfully applied to two dimensional 
and three dimensional problems with uniform mean flows. 

For swirling mean flows, centrifugal and Coriolis effects produce force imbalance and, 
as a result, the acoustic, vortical and entropic modes are coupled [6, 7]. Moreover, the 
eigenvalue problem resulting from the normal mode analysis is not self-adjoint and, as a 
result, the modes may not form a complete set. Another consequence of the force imbalance 
is the development of a singular critical layer for the vortical modes, where it is very 
difficult to obtain accurate numerical solutions. To overcome this difficulty, the normal 
mode representation is limited to the acoustic pressure-dominated modes which propagate 
with a phase velocity equal to or larger than the local speed of sound in a reference frame 
moving with the local mean velocity [8]. In order to obtain a complete representation of the 
disturbance field, an initial-value-analysis is used to determine the vorticity and entropy- 
dominated fields which propagate with phase velocities close to the mean flow velocity. An 
iterative scheme which enforces conservation of mass momentum and energy is then used to 
ensure that the superposition of the nearly sonic and nearly convected fields exactly satisfy 
the Euler equations. 

In section 2, we present a normal mode analysis to represent the acoustic modes in 
an annular duct and a criterion based on group velocity is used to satisfy the causality 
condition. In section 3, we present an iterative scheme for the inflow/outflow conditions 
and implement this scheme in an Euler code. The boundary conditions are then tested on 
the propagation of harmonic disturbances from an upstream imposed source in an annular 
duct with swirling mean flow. In section 4, an expression for the sound power in swirling 
flows is derived which may be used as a conservation law and to study the effect of mean 
swirl on the sound power. 
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2 Normal Mode Analysis in Nonuniform Swirling Flows 


This section presents the normal mode analysis of the linearized Euler equations for nonuni- 
form mean flows. The outgoing acoustic modes are used in the next section to formulate 
and implement nonreflecting boundary conditions. 

2.1 Introduction 

Sound propagating in a duct is usually expressed in terms of the duct normal modes. A 
major advantage is that at a given frequency, there are only a finite number of propa- 
gating modes. In aeroacoustic computations, the modes can also be used to derive exact 
expressions for the nonreflecting inlet /outlet boundary conditions. 

As mentioned earlier, the mean flow velocity of fan engines is characterized by significant 
swirl and radial variations in both the axial and swirl components. Moreover, the work done 
by the upstream fan blades produces radial variations in the flow enthalpy, and losses in 
the boundary layers near the hub and tip walls produce significant changes in the entropy. 
Representative results from a Reynolds averaged Navier Stokes calculation are shown in 
Figure (1) for the circumferentially averaged total enthalpy variations downstream of the 
fan. Note the large rise in the total enthalpy in the tip region. This rise is associated with 
the tip vortex which is shed downstream of the fan. 

Current fan noise schemes compute the wave modes in a constant area annular duct by 
assuming the mean flow is axial and uniform [9] or that the annulus is narrow [10, 11]. The 
normal mode analysis can then be used to determine stator vane counts which cut-off blade 
passing frequency tone noise [12]. 

In order to incorporate nonuniform mean flow effects, a normal mode analysis of such 
flows must be carried out. Recent studies, to account for swirl, have assumed simple analytic 
models for the mean flow swirl which consist of combinations of free vortex and rigid body 
swirl [13, 14, 15]. These treatments assume uniform total enthalpy and entropy. The 
results show that the centrifugal and coriolis effects create force imbalances which couple 
the acoustic, vortical and entropic flow disturbances. However for the high frequencies 
relevant to fan tone noise, the coupling is weak ([14], [8]), i.e. the acoustic modes have 
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small vorticity content and the vortical modes have small pressure content. 

Another important application of the normal mode analysis is the formulation of non- 
reflecting boundary conditions. This, however, requires an accurate solution of the wave 
modes in a duct [3]. To examine the sensitivity of aeroacoustic computations to the accu- 
racy of the acoustic modes used in the formulation of the nonreflecting boundary conditions, 
we consider the case of , a one-dimensional incident acoustic wave whose pressure is of the 
form with an imposed outlet boundary condition error in the mode wavenumber, 

A k. The solution is given by, 

u(z, t ) = [e ikx + Sm ^ri e ikL (iAkL + 0(AkL) 2 )]e~ iut . (1) 

sin (kL) v 7 

Note that the error in the numerical solution is linearly proportional to the error in the 
downstream boundary condition, A k. Moreover, the error becomes very significant if the 
duct length is an integer multiple of A/2. 

Figure (2) shows the numerical solution and error in a domain of length .75 which is 
discretized with 33 points. The real part of the numerical solution is shown in the bottom 
plot of Figure (2) for various boundary condition errors, A k = 0.0,0.025,0.05,0.1, on the 
wavenumber, k = 2i r. The relative error of the solution is shown in the top plot in Figure 
(2). Note the relative error of the numerical solution is directly proportional to the error 
at the exit boundary. For example, as the exit boundary condition error increases from 
.025 to .05 the relative error of the numerical solution increases from nearly .11 to .22. 
This degradation of the solution is visible in the solution where the exact solution is the 
solid black line. As the boundary condition error increases, modifications to both the 
phase and amplitude of the numerical solution are visible across the entire domain. The 
apparently small error right at the exit boundary, as is indicated by the small difference 
between the exact solution and the numerical solution at the end of the domain, is because 
we have plotted the real part of the solution and A k is relatively small. This simple result 
highlights the importance of accurately specifying the exit boundary condition for general 
axisymmetric mean flow to obtain an accurate numerical solution. 

In what follows, we (i) generalize the duct mode calculation for axisymmetric mean 
flow to allow arbitrary radial distributions, (ii) examine the effects of arbritrary swirl and 
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enthalpy on the acoustic modes in a duct and (iii) compare this with the uniform flow, 
narrow annulus and uniform enthalpy models which have been used in earlier work. Many 
of the Results for the normal modes with arbitrary radial distributions are also in Ali et al 
[16]- 

2.2 Mathematical Formulation 


For a non-viscous, non-heat-conducting perfect gas with constant specific heats, the gov- 
erning equations are the Euler Equations, 


2-U + U-VU = yVp 

= O' 

Dt 


( 2 ) 


The flow variables are decomposed in terms of the sum between their steady mean values 
and their unsteady perturbations, 


U(x,t) = Uo(x) + u(x,t), 
p(x,t ) = p 0 (x) + p'(x,t) 

p(x,t) = po{x) + p\x, t), 
s(x,t) = 5 0 (x) + s'(x,t), 


(3) 

(4) 

(5) 

( 6 ) 


where x represents any coordinate system, Uq, po, Po-i So are the steady mean velocity, 
pressure, density and entropy, respectively, and u, p 1 , s' are the corresponding un- 
steady perturbation quantities. The unsteady quantities are assumed to be small such that 
|u(T, 1)|, |p'|, |p'|, l^'l < Uo{x), po, po, So- We non-dimensionalize with respect to the mean 
radius of the duct, r m , and the mean density, p m , and speed of sound, c m , at the mean 
radius. 

Substituting into (2), the nondimensional Euler equationsg reduce to the following cou- 
pled system of linear equations. 

^{u - s'/(2c p )U 0 \ + [(u-s'/(2c p )Uo) • V]f/ 0 = -^V(-^r), (7) 

L>t Po p l 0 n 


—' J £ 1 -s'H2c p )}=-±-V-(poS), 


Dt p 0 cl 


Po 


( 8 ) 
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(9) 


^s' + u-VS o = 0, 


where is the mean flow material derivative 


Tt s l + (Vv ' 


( 10 ) 


The impermeablity condition is imposed at the hub and tip radii. 


u r = 0. 


(id 


2.3 Mean Flow 


The mean flow is assumed axisymmetric and, away from the upstream and downstream 
blade rows, the axial mean flow gradients are treated as much smaller than the radial 
gradients. In such a flow region, mass conservation implies that the radial velocity is small 
and the nondimensional mean flow velocity can be expressed in the form, 

Uo(x) = M x (r)e x + M s (r)eg , (12) 


where e x and eg represent unit vectors in the axial and circumferential directions, respec- 
tively and M x (r) = U x /c m and M s (r ) = U s /c m . The mean flow is, in general, vortical with 
vorticity given by 


t „ ft 1 d(rM s ) dM x 

(o = V x P 0 = -j e x -j— eg. 

r dr dr 


(13) 


The stagnation enthalpy, entropy, velocity, and vorticity are related by Crocco’s equation 


VH = TV So + U 0 x Co- 


(14) 


In order to study the effect of the mean flow on the eigenmodes in an annular duct with 
swirling flow and to test the robustness of the current scheme, two mean flows are used: 

a ) Theoretical Model represented by a combination of rigid body and free vortex swirl, 

M s = M u + M r (15) 

where M,^ a r and Mr a 1/r. Assuming uniform enthalpy from hub to tip and an 
isentropic flow, the axial component of velocity takes the form, 

Ml = M 2 xm - 2 [ft 2 (r 2 - 1) + 2ftrin(r)], (16) 
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where M xm is the axial Mach number at the mean radius of the duct, 
b ) Mean flow data from a fan design. 

Mean flow data representative of a fan design [17] is compared to the first theoretical 
mean flow model. Figure (3-a) shows the total temperature which is proportional to the 
stagnation enthalpy across the annulus. For this fan design, the stagnation enthalpy is 
not constant but it varies 8% from hub to tip while the theoretical model gives constant 
stagnation temperature. Figure (3-b) shows the swirling component of the mean flow data 
along with the best fit of equation 15. Using the least square method, Mn m = 0.2108 and 
Mp m = 0.2498. In figure (3-c), the axial component of the mean flow data is compared to 
the values obtained from equation (16). The figure shows a maximum difference of 15% 
between the two curves. 


2.4 Normal Mode Expansion 


Normal mode analysis is used to obtain the spectrum of propagating and evanescent 
acoustic-vorticity modes for the swirling flow equations. The following Fourier expansion 


is assumed 

{u x , u r , ue, p, s'}(z, r, 6; t) 


roc oo 

/ E 

^ 00 m=— oo n=l 

T mn {r),p mn {r),s mn {r)} x (17) 


where m, and n are integer modal numbers characterizing the circumferential and radial 
eigenmodes, respectively. Since the equations are linear, each Fourier component can be 
considered /eparately. The coupled system of equations (7,8) and (9) can now be repre- 
sented in terms of normal modes as 

Mr . dM x , . , . 1 . 

A mn(Xmn ~ 7J s mn ) H ~j ( - 1 ^ran) = ~^mn(~Pmn)i (18) 

2 c p dr po 


iA-mnRmn 4 " 


2 c p r' 


2 M, 


~ Smn. 


S T 4- 

L mn i 


1 dpmn _ J7 
Po dr (iPoPo) 


Pmn — 0, 




1 d 


m 


A mn{T mn — — S mn ) + -—(rM s )( iRmn ) + n ^ mn ~ 


r dr 


Por 


(19) 

( 20 ) 
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— S mn /(2 c p )) + —T mn + -j^—[po r ]( — iR m n) + iRmn) — k mn X mn , (21) 


po4 


dS n 

dr— iff 


'-'TTin — ^ 

where the convective eigenvalue, A mn , is defined by the expression, 

mM s 

Amn = — ^ 4" k mn JvI x 4" • 


( 22 ) 


(23) 


The boundary condition at the hub and tip radii is 


Rmn(r) = 0. (24) 

Note that unlike the uniform axial mean flow limit this is not a Sturm-Liouville eigenvalue 
problem. Therefore there is no proof of completeness or orthogonality of the eigenfunctions. 


2.5 Numerical Methods 

A combination of a pseudo-spectral method and a shooting method is used in order to 
obtain the eigenmodes and eliminate any spurious modes. In the pseudo-spectral method, 
the eigenfunctions are represented by a series of Chebyshev polynomials which form a 
complete orthogonal set [14]. 

For an annular duct, the physical domain < r < r t is mapped to -1 < x < 1. 
Choosing the grid points as x\ = cos[?r(/ - l)/(iV - 1)] for 1 < / < N, the eigenfunctions 
can be represented as, 

/, = £>,«*[*(„- 1)(/- 1)/(JV- 1)], (25) 

V - 1 

where a p are coefficients of the spectral expansion. Derivatives then are given by // = 
and /" = where and are known coefficient matrices [18]. 

The problem is reduced to solving an algebraic eigensystem of the form, 

[A]x = k mn [B]x (26) 

where the matrices [A] and [5] result from equations (18-21) after discretization and x is a 
vector of the eigenfunctions X mn , p mn , Rmn, and T mn . The algebraic eigensystem is then 
solved by a shifted QR algorithm [19]. 
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This method has exponential convergence and can resolve thin layers of steep changes. 
The disadvantage of this method is that the matrices [A] and [B] are usually not sparse 
which may lead to difficulties in solving the eigensystem as the number of grid points 
increases. It also has the disadvantage, like the finite difference method, that this system 
of 4N equations has 4 N solutions. Since the number of physical solutions in a problem is 
independent of the numerical grid, the vast majority of these solutions are spurious. 


An alternative method to solve the eigenvalue problem is the shooting method. In this 
method, equations (18-21) are cast in the following form, 


Zr Y = F(Y) - 


(27) 


In the present problem, Y — {p m n,—iRmn} T •> and equation (27) becomes, 

7 Pmn 


dpmn - o n (-iA R - + 2— T ) + -fe 

dr 2 c p r r (7 p 0 ) 


(28) 


d{ iRmn) . /Pmn 

— 2 


dr 


Poc 0 


o m illPo) 

4 -k x 1 <t 

) ± mn 1 ^mn^mn i 

c p r p 0 r 


(l iRmn )• (29) 


The axial and tangential velocity components T mn , X mn are functions of p mn and R mn , 

(30) 


_ 1 r/ • r> d ci / m , M s 

Tmn ~~ T ~ 7 ~V M s ) Pmn\ "b ^ 

Amn r dr p 0 r 2 c p 


A mn — 


A mn df 

and s mn is determined from (22). 


1 \ dMx (jR \-hl Hr, +k —S ] + — S 

Po 2c p 2c p 


(31) 


The initial conditions at are given and assuming an initial guess for fc mn , the previous 
equations are integrated using a Runge-Kutta scheme from r ^ to r t . The boundary condition 
at r t is checked. If this boundary condition is not satisfied, another value of k mn is assumed 
until the boundary condition at r t is satisfied. Newton method or secant method is used 
for the iteration process of k mn . The main advantage of the shooting method over the 
finite difference and spectral method is that it does not produce spurious modes. However, 
it needs an initial guess for each eigenvalue which makes it impossible to cover all the 
eigenvalues in the complex plane. Also, it may not directly yield solutions with discontinuity 
or very steep gradients. In the present work, we use a combination of the pseudo-spectral 
method and a shooting method. The pseudo-spectral method provides the initial guess for 
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the shooting method which eliminates the spurious modes and improves the accuracy of 
the solution. 

Figure (4) shows a comparison of the eigenmode spectra obtained by the pseudo-spectral 
method and that obtained by the shooting method for a case of m = — 1 and u — 16 using 
the isentropic mean flow data in figure (3). The vertical axis corresponds to the imaginary 
part of the axial wavenumber and the horizontal axis corresponds to the real part. Those 
modes which lie on the horizontal line corresponding to zero imaginary part represent un- 
steady disturbances which propagate downstream without decay. The modes with positive 
imaginary part correspond to evanescent modes. A typical eigenmode spectrum for swirling 
flows consists of two families of solutions. The first family represent the pressure-dominated 
modes, herein referred to as the acoustic modes, and consists of propagating modes and 
evanescent modes. This family appears on the left side of figure (4). The other family 
consists of vortical modes and appears on the right side of figure (4). A calculation of the 
group velocity shows that the acoustic modes to the left of the evanescent modes prop- 
agate upstream and the acoustic modes to the right of the evanescent modes propagate 
downstream. The group velocity was calculated numerically since there is no closed form 
expression for it. In figure (4), two of the acoustic propagating modes obtained by the 
pseudo-spectral method are spurious; one of them appears to propagate upstream and the 
other appears to propagate downstream. However, using these eigenvalues as initial guesses 
for the shooting method results in a single acoustic mode. Some of the evanescent modes of 
the pseudo-spectral method also are inaccurate. It is difficult to assert using the shooting 
method whether the nearly convected modes are spurious or not. This is due to the singular 
behavior of the eigensolutions at certain values of the radius. However, due to the singu- 
lar behavior of the eigensolutions in the critical layer an initial value approach is derived 
in section 3 to completely represent the non-acoustic disturbances for the inflow/outflow 
boundary conditions. In section 2.6, we compare our eigenvalue algorithm against the nar- 
row annulus limit and the uniform enthalpy limit to assess the range of validity of these 
cases. 
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2.6 Narrow Annulus Limit 


In the narrow annulus limit, analytic expressions for the dispersion relation and radial 
eigenfunctions can be derived. In order to examine the accuracy of the present model 
and assess the limitations of the narrow annulus model, we briefly present the dispersion 
relation and eigenfunctions of the acoustic modes in the narrow annulus limit and compare 
the number of propagating modes and the radial eigenfunctions obtained analytically with 
those determined from computation of the general problem. 

Since the annulus is assumed narrow, i.e. r^jvt << 1, we can expand the radius, 
mean velocity, density, pressure and speed of sound about the mean radius in terms of a 
perturbation expansion of the form, 


r — 1 + ef + • • • 

U 0 — M x (l)e x + M$(l)ee + 0(e) (32) 

(po-i Poi Co) — ( i i Poi 1 ) 0(e) 


where |c| << 1. As a result, the leading order mean flow quantities are uniform. In this 
limit, the linearized Euler equations become, 

!§£ - VV = o (33) 


where ^ = JL + + M s -^. To satisfy impermeability at the hub and tip, we impose 

dp' n 

— = 0, r = r h ,r t . (34) 

Solutions to (33, 34) exist of the form, 


m — oo n — oo 

P'= E E ™s[lmn(r - r h )]e i ( me+k ™ 1 -^ (35) 

m= — oo n = 1 


where 


-&M X ± - /? 2 (m 2 + 7g, w ) 

^2 


(36) 


with J; = uj-mM s ,f3 = 1 -Ml and 7 mn = ( r J!i rh y Equation (36) indicates that the number 
of propagating modes depends upon the reduced frequency of the disturbance. The swirl 
component of velocity acts as a Doppler effect to increase the frequency if the disturbance is 
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rotating in a direction counter to the swirl component and to decrease the frequency if the 
disturbance rotates in the same direction as the swirl. Moreover, the number of propagating 
modes grows with frequency. The analytical solution also indicates that if a radial mode 
is propagating at a given frequency, it will propagate at higher frequencies. As a result, 
one can determine all the radial modes at all frequencies by computing the radial modes at 
the maximum frequency of interest, u = Cj max . In what follows, we examine solutions for 
both narrow the annullus and the general case to examine the extent to which the narrow 
annulus approximation is adequate for real three dimensional geometries. 

2.7 RESULTS and DISCUSSION 

In this section, we present the normal mode solutions for a variety of mean flows and hub- 
tip ratios. The eigenfunctions, p mn , are normalized numerically such that the norm of each 
mode satisfies, 

/ [r\\Pmn\\ 2 ]dr = 1. (37) 

Jr h 

2.7.1 Effect of Hub-Tip Ratio and Frequency on Duct Modes 

We examine the effect of hub-tip ratio and reduced frequency on the duct modes. For 
simplicity, we consider a mean flow with free vortex swirl and examine two cases. The 
problem is characterized by the parameters u r = ur m /cQ , r* = r^/r* and f = Y 

In case one, we consider a large spinning mode order, m — —20 with u r = 30, r* = .4 
and f = .21. Figure (5), shows there are 6 downstream propagating radial eigenfunctions. 
Each radial mode is characterized by a different number of zero crossings. We note that the 
modes are quite different from the sinusoidal solution obtained by the analytic result. For 
example, the amplitude of the numerically determined eigenfunctions are strongly weighted 
in the hub region and only become sinusoidal in the tip region. This is not surprising since 
the hub-tip ratio for this case is outside the range where the narrow annulus solution would 
be expected to be valid. 

In case two, we consider a smaller spinning mode order, m — — 1, with r* = .4 and 
T - .21 and compare the solutions for two different reduced frequencies, uj r — 20,40. In 
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figure (6), we plot the sixth radial mode for the two reduced frequencies taken in case two. 
Note that, as in the narrow annulus limit, the radial eigenfunctions are independent of the 
frequency. 

2.7.2 Effect of Mean Flow Model 

Figure 7 compares the eigenmodes obtained using the mean flow data and the theoretical 
mean flow model of equations (15) and (16) for m = -12 and u = 24.2. The number of 
acoustic propagating modes is six modes in each direction for the two swirling flow models. 
For the case of no swirl, the number of acoustic propagating modes does not change. The 
effect of the mean flow model on the downstream acoustic modes is small while the effect 
is more significant on the upstream acoustic modes. The strongest effect of the mean flow 
model is on the convected modes. 

The effect of the mean flow model on the eigenfunctions of the six upstream acoustic 
propagating modes is shown in figure (8). The mean flow model has a significant effect on 
the eigenfunction radial variation. As the order of the mode increases this effect decreases. 

The effect of the mean flow model on the eigenfunctions of the six downstream acoustic 
propagating modes is shown in figure (9). The effect of the mean flow model on these 
eigenfunctions is significant only for the first two modes. This result is not surprising 
since it is expected that the longer wavelength radial disturbances will be more influenced 
by the mean flow gradients which vary on the mean radius lengthscale than on the short 
wavelength higher order radial modes. 

2.8 Summary 

A numerical method was developed to compute the normal modes in an axisymmetric mean 
flow with arbitrary swirl. Two sets of propagating modes were found. One set is pressure- 
dominated and propagates with a phase velocity, in a reference frame moving with the 
mean velocity, greater than or equal to the local speed of sound. The other set propagates 
with a phase velocity close to the mean flow velocity. Numerical results were presented for 
a variety of mean flows and compared with the narrow annulus and uniform enthalpy limits 


NAS A/CR— 2003-2 12311 


13 



to examine the importance of centrifugal effects and enthalpy variations on the acoustic 
modes. The results lead to the following conclusions, (i) The number of propagating modes 
is strongly dependent upon whether the modes are co-rotating or counter-rotating to the 
mean swirl due to the doppler shift in frequency, (ii) Mean flow swirl refraction effects are 
more important for the upstream-going modes than the downstream-going modes resulting 
in larger deviations in the eigenvalues of the upstream-going modes. This is because the 
upstream modes propagate with lower phase speeds, while the downstream modes propagate 
with higher phase speeds, and are thus less affected by the mean flow swirl velocity, (iii) The 
narrow annulus limit misses centrifugal effects which are important for moderate hub-tip 
ratios and as a result misses the behavior of the eigenfunctions in the hub region, (iv) The 
uniform enthalphy limit is a reasonable approximation for the higher order radial modes 
but is inaccurate for predicting the behavior of the lower order radial modes. 

The generalized eigenvalue scheme presented is a necessary module of any computational 
acoustics scheme which employs nonreflecting boundary conditions. The coupling of this 
analysis to a numerical scheme is presented in section 3. 

3 Boundary Conditions for Aerodynamic and Aeroacoustic 
Computations 

3.1 Introduction 

In this section, we derive time-harmonic inflow/outflow conditions for unsteady aerody- 
namic/aeroacoustic computations based on the linearized Euler equations. We first exam- 
ine the results of a uniform mean flow which was treated and applied to various three- 
dimensional wake and acoustic wave excitations in [3]. In this case, at the inlet of the 
computational domain, the inflow conditions imposed an incident acoustic or vortical wave 
and applied a nonreflecting condition for the radiated pressure. At the outlet boundary, 
we specified two conditions: a nonreflecting boundary condition for the pressure and a 
pure convection condition for the vortical part of the velocity. The exact outlet conditions 
result from the property that in a uniform mean flow the acoustic, entropic and vortical 
disturbances can be split into distinct modes [5]. 
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For mean flows with swirl, it is not possible to split the disturbances into distinct modes 
because the local Coriolis and centrifugal forces couple the various modes. Thus, as was 
shown in the previous section, there are pressure-dominated modes which mainly represent 
the acoustic part of the unsteady flow and vorticity-dominated modes propagating at nearly 
convected velocities. Moreover, there is an infinite number of vorticity-dominated modes 
which cluster at each side of the critical layer [7]. 

Montgomery et al. [4] attempted to circumvent this difficulty by approximating the 
vorticity-dominated disturbances as pressure-free convected disturbances. However, this 
representation does not satisfy the radial momentum equation in nonuniform flow. Gol- 
ubev &: Atassi [8, 20], used an initial- value analysis to represent the vortically-dominated 
disturbances and have combined this analysis with the normal mode analysis to develop 
accurate boundary conditions for high frequency, rigid body and free vortex swirling flows. 
Their results, along with Elhadidi et. al [21], show that the stability and growth of the 
vorticity-dominated modes depend on the mean flow swirl distribution and affects the up- 
wash velocity along a blade surface. 

In the present work, we extend the analysis of [8] to general swirling flows and present an 
algorithm for the boundary conditions which accounts for the coupling between the acoustic 
and vortical modes. We implement the boundary conditions in a linearized Euler solver 
and, as an application, we test the boundary conditions by computing the propagation of 
acoustic and vortical disturbances in an annular duct. Sensitivity of the solution to inflow 
conditions is examined and, in particular, the effect of representing the inflow excitation as 
a convected disturbance is studied. 

In section 3.2, we derive nonreflecting inflow/outflow boundary conditions for swirling 
flows. The Euler scheme and the propagation of acoustic and vortical disturbances in an 
annular duct with nonuniform mean flow is presented in section 3.3. In section 3.4, we 
compute the propagation of vortical and acoustic waves in a variety of mean flows and 
examine the importance of the inflow wake specification on the evolution of the wake. 
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3.2 Nonreflecting Boundary Conditions 


Based on the normal mode analysis in section 2, the pressure-dominated modes can be 
expressed in terms of the eigenfunction expansion, 


771=00 n=oo 


pi = 53 

77l=— OO 
771 = 00 

E “i-pUf)'"*"'*"' - "' 1 . 

71=1 
71 = 00 

(38) 

p7 = 


(39) 


ra= — oo 7i=l 


where the + (-) indicates modes propagating with a group velocity downstream (upstream). 
Thus, pf corresponds to the sonic part of the pressure whose energy propagates downstream 
and pj represents the acoustic pressure whose energy propagates upstream. 

In order to avoid small terms in (38) which may cause numerical difficulties, we use a fi- 
nite number of propagating and slowly decaying modes to represent the pressure. The series 
expansion can be truncated to a finite number of terms with negligible error. Discretizing 
and truncating the expansion so that there are N r radial modes and Nq circumferential 
modes we obtain, 

m=Ne/2-l n =Nr 

pf ]k (x (h 6 jk ,r jk ,t) = £ £ a+ n pi n (r jk )e l ^o + rne ]k -.t)^ (40) 

m--N e / 2 n=l 
m= Ng/2-1 n=N r 

P7 jk (XO,0] k ,r jk ,t) = £ E (41) 

m = -N e /2 n= 1 

Note that p± k can be cast as a vector with values at each of the inlet plane points (Ojk, T'jk) 
in the form, 

{pi,}=(JO ± )e~“‘. (42) 

where are N r Ng x N r N$ square matrices containing the axial, circumferential and 
radial eignfunctions for each of the N r radial modes and N$ cicumferential modes. The 
amplitude vectors containing the amplitudes, a^ n , of the N r Ng modes is denoted by A ± . 
The pressure is evaluated at N r No points in the inlet/outlet plane. The indices j,k refer 
to the grid points in the axial plane, x = xq , with j = 1 ,2,..., Nq and k — 1 , 2, ..., N r . 


NASA/CR— 2003-2123 1 1 


16 



3.2.1 Inflow Boundary Condition 


The unsteady pressure, velocity and entropy can be expressed as the sum of an outgoing 
nearly sonic field and a known incident field which is some combination of acoustic and 
vortical waves, 

P - P7 +pI + pL (43) 

u = uj + u{ + u l c , (44) 

s' = sj + s l s + s' c . (45) 

We assume that the known incident field, p 1 = p[ + pl,u ! = u[ + u^s 1 = si + sf., is 

obtained via experiment or computation. Given this information, we propose the following 
boundary condition at the inlet boundary. 


Nonreflecting Boundary Condition Algorithm 

1. The sonic incident field, p[ is obtained by projecting p l onto the incoming normal duct 
modes. With p[ determined u [ and si are determined from the normal mode analysis given 
in section 2. 

2. The nearly convected incident field is then given by = u 1 — u^pl = p 1 —pl, si = s 1 —si. 

3. The inlet condition for the upstream-going acoustic waves, pj, is expressed 

P7 = (p- P 1 ) = M 7o M (xo\* z)(p - P%o+*x)- (46) 

Given pj, we can determine from the normal mode analysis, given in section 2, uj and sj. 

3.2.2 Outflow Boundary Condition 

At the outflow boundary, the unknown solution contains a superposition of outgoing acous- 
tic and vortically-dominated disturbances, 

P = pf + Pc, 

u = uf + u c (47) 

s' - sf + s c . 
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In contrast to the inlet boundary, the outlet boundary has two unknowns: the downstream- 
going acoustic waves and the downstream-going nearly convected disturbances. As a result, 
two conditions are needed to specify the outlet boundary condition. One for the outgoing 
acoustic waves and one for the outgoing nearly convected waves. 

Nonreflecting Boundary Condition Algorithm 

1. The exit condition, x = xjy x , for the outgoing acoustic pressure waves is expressed by 


~+ _ 


= (p-p^) = M+ M+- . 




(48) 


Initially, we set p° ld = 0 implying pf = p. In later iterations, we express pf = p — p° ld - 
Given p+ we can determine from the normal mode relations in section 2 ut and sf . 

2. The nearly convected part of the solution is given by u c — u - uf, s c — s' - sf . The new 
convected pressure field, p ™ ew , is determined from using the mass and momentum equations 
(7,8). A condition of the form + u c -VU = 0 is applied on the convected velocity field. 

3. Check p™ ew — p° c ld . If difference is greater than 1 x 10 -2 , set p° ld — p™ ew . 


3.2.3 The High Frequency, Homentropic Limit 

The exit boundary condition initially assumes that the nearly convected pressure field, p c , is 
small and then iterates in time to reach a converged solution. In order to justify this initial 
guess, we consider the high frequency and homentropic, limit. In gust response applications, 
the reduced frequency is large and the radial entropy gradients are small except near the 
tip of the annulus where the tip vortex lies. In this limit, significant simplifications to the 
exit condition can be attained. 

The vorticity-dominated disturbances are periodic over the upstream blade row passage 
and move with a phase velocity close to the mean flow velocity. Without loss of generality, 
these disturbances can be expanded in the form [8], 

m= oo 

(u c ,p c ,s c )= ^2 {^rn(x,r),p c m (x,r),s c m (x,r)}e ,4>m , (49) 

m — — oo 

where o m — a m x + md — ut is defined such that 

Do , 

~ (50) 
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and m is an integer index which for tone noise is a multiple of the number of blades upstream 
and A m = ( A mx , A m g, A mr ). Equation (50) implies that 


— 


(u> — M$ m/r ) 

M~ x 


(51) 


In the high frequency limit, the gust disturbance scales with the blade passage lengthscale, 
l g = (27rr m )/(m) and the mean flow variations scale with the mean radius. Typically, the 
number of fan blades is large implying that for tonal applications l g /r m is small, 0(c), and 
the reduced frequency, is large, 0(1 /e) where e << 1. We further assume that a m is large 
from (51). In this case, the nearly convected disturbances propagate as slowly modulated 
wavetrains with rapidly varying phase. As a result, the leading order contributions to the 
Euler equations become, 


xr d a m r 

( ) A mr T A m Q T Amx — 0? 

m dr m 


xr dA mx dM : 
M x — — + 


dx 


dr 


x A - 


P 0 




ox r dr p 0 r 


M, 


dA 


mr 0 . Id 

- 2 —A m e = -a m (ixp c m ) 


(52) 

(53) 

(54) 

8x " r “ m</ pQ dr~" lK ~ rm ' (55) 

where although a m , m are large (a m p c m ) and (mp c m ) are 0(1). This implies that in the high 
frequency limit, p c m is 0(c). Note that even in the case where only a m or m is large p c m is 
0(c). As a result, the pressure field associated with the vortically-dominated disturbances 
is small and the pressure is to leading order determined from the sonic pressure, p s . Thus, 
in the high frequency limit, the corrections to the pressure field resulting from the iterative 
process are small and the vortical and acoustic disturbances are uncoupled. This result 
suggests that the exit boundary condition can be satisfied to leading order with only one 
interation. 


3.3 Wake Excitation 

In th is section, we calculate the wake excitation from the rotor-locked flow field. In previous 
work [4], experimental data was taken for only two components of the velocity and so at 
the inflow the wake was modeled as a convected disturbance with zero radial velocity and 
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pressure. Such an inflow specification does not satisfy the linearized radial momentum 
equation. Moreover, the second blade passing Fourier component of the rotor locked field 
shown in figure (15) shows that the velocity field is fully three dimensional with a significant 
radial velocity component. In order to uniquely specify the wake four of the five (three of 
the four for homentropic flow) flow variables must be measured or computed via a CFD 
calculation of the fan to uniquely specify the inflow source, u ! ,p ! and s 1 . The unknown 
component can then be determined from the continuity equation. 

To determine the inflow excitation, the rotor-locked flow velocity far downstream of the 
upstream blade row is assumed periodic over the passage of the upstream blade row. As a 
result, the flow variables can be expanded in a Fourier series, 

771 = 00 

E <e - tmN » 6 

m=— oo 

771 = 00 

P r = £ P I me- imNb e (56) 

77l = — OO 

771 = 00 

*'= £ s I me - imNb9 


where the Fourier coefficients are given by 


27 r Je 

s' m = ^ [“*% s'e^'de'. 
27t Je 


(57) 


The flow variables u / ,p / ,s / can be determined from Reynolds averaged Navier-Stokes cal- 
culation or experiment. For fan applications, the entropy variation is typically small except 
in the tip regions. As a result, the inflow plane may be assumed homentropic and only- 
three components of the velocity need to be computed or measured. 


3.4 Numerical Formulation for the Interior Domain 

In this section, we assume the flow to be governed by the linearized Euler equations (7, 
8. 9) and we briefly present the interior Euler scheme for computing the propagation of 
acoustic and vortical disturbances in a duct. This model problem enables us to examine 
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the convergence and accuracy of the numerical solution and assess the accuracy of the 
nonreflecting boundary conditions. The characteristic lengthscale, L = r m , is taken to 
be the mean radius of the annulus. Nondimensionalizing, assuming solutions of the form 
V = U(x, t)e~ iwt , the linearized Euler equations can be expressed 

Cl! = -f . (58) 

where U = [/3, u x , ug, u r , s] T , C = — I(w*>) + A x ^ + + C r ^ + The variable coef- 

ficients A x - D depend upon the mean flow and the mean flow gradients. Equation (58) 
is discretized with a second order accurate Lax Wendroff scheme. Note that although the 
problem is time-harmonic, the time derivative remains and represents whatever transient 
time-dependence exists from the initial condition. The term ^ or residual is used to denote 
that the numerical solution iterates in time to a steady state solution. Hence, when ^ -* 0 
a steady state solution is obtained consisting of a single frequency, u. The impermeability 
condition at the hub, r = rh/r t , and the tip, r = 1, of the annulus imply 

u r (x, 6, rh/rt) = u T (x, 6, 1) = 0. (59) 

To relate this model problem to harmonic rotor-stator interaction whereby an upstream 
rotor consisting of Nb blades interacts with a stator consisting of N v structural vanes, the 
solution is assumed quasi-periodic in 9. In this case, the computational domain extends from 
[0,2 ir/N v \. Quasi-periodicity implies U(x,0, r) = U(x,2ir/N v , r)e l(T where a = 2wnNb/N v 
with n being the harmonic of the upstream excitation. To complete the formulation of this 
initial-boundary-value problem, nonreflecting boundary conditions, as described in section 
3.3, are imposed at the inflow and outflow boundaries of the computational domain. 

3.5 Numerical Results 

The Euler scheme is applied to compute the propagation of acoustic and vortical waves for 
a variety of mean flows ranging from uniform mean flow to simple swirling flows consisting 
of some combination of free vortex and rigid body swirl. We also consider a variety of inflow 
excitations ranging from an acoustic excitation to a vortical excitation to a combination of 
the two. Since there is no scattering body, no reflected waves should reach the upstream 
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boundary. The upstream nonreflecting boundary condition is then, 


U(0,e,r)=U I (0,e,r). 


( 60 ) 


3.5.1 Uniform Mean Flow 

The u nif orm flow test case is motivated by the problem proposed by the CAA workshop 
[22]. An incident wave consisting of both a vortical transverse disturbance and a pressure 
disturbance is imposed for various hub- tip ratios and radial phase variations. The incident 
disturbance is, 

U 1 = A(r)e im 6 + B(r)e i[Nbl>+2lTq{ ^ )] (61) 

where A(r) is the amplitude vector of the acoustic disturbance and B(r) is the amplitude 
vector of the vortical disturbance. The amplitude of the incident disturbance is given by 

A t = [*(w- M g k+ n )p+ n ,ik+^+ nt (im/r)p+ n , ^Pmn], 

B t = k[ 0, ~ MriVfe ,l,0] (62) 

ur 

where n is used to indicate the size of the vortical disturbance relative to the acoustic 
disturbance. The parameters q and r^/r* give the degree of radial phase variation in the 
vortical gust. The outflow boundary condition (46) is used for the remainder of the infinite 
domain. Note that in the uniform flow limit one only needs to express a condition for the 
pressure [3] since the acoustic and vortical fields are decoupled. For simplicity, a uniform 
grid is used to study the convergence of the fully three-dimensional case. The range of 
parameters studied in the five cases are shown in table 1. To obtain a measure of the 
convergence of the scheme and its evolution to a steady state the average relative error, 

p v r; - Pill , m real 

E - 1 “iMr /(iV) (63) 

is computed at each timestep where N is the total number of points in the domain and i 
is the grid index. Upstream, we impose a disturbance, [/^(O,#). From the normal mode 
analysis the ‘exact’ solution, {7 e , can be obtained. We consider five cases to study how well 
the numerical solution converges for problems with varying degrees of three-dimensionality. 
The radial phase of the gust is the parameter which governs the strength of the radial 
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Parameters 

Case I 

Case II 

Case III 

Case IV 

Case V 

uL/c 

2tt 

2 ir 

2r 

27 r 

2tt 

a 

47t/3. 

4n/3 

47t/3 

47t/3 

4tt/3 

M x 

.5 

.5 

.5 

.5 

.5 

rh/rt 

.65 

.65 

.65 

.65 

.65 

Q 

1/4 

1/2 

i 

i 

2 

K 

.1 

.1 

.1 

i 

1 


Table 1: Selected parameters for five cases of wave propagation in a uniform mean flow. 

variation in each case. The parameters, q and k are varied, while the other parameters 
are held constant in cases I-V as shown in table 2. In each case, the variables (p,ug,u r ) 
are plotted against x, 6 and r respectively and the global relative error of the density is 
plotted versus time iteration. In Figures 10, 11 and 12, the solution and the relative error 
are plotted for cases I, II and III. The grid in all the cases considered is 80 x 25 x 20 in 
x, 6 and r, respectively. The dashed lines correspond to the exact solution and the solid 
lines correspond to the numerical solution. The agreement is quite good and the relative 
error of the density is .1 percent when the solution reaches a steady state. In this case, 
q = 1/4, 1/2, 1 with k = .1 and the increasing variation in r is most observable in the plot 
of ug versus r. 

In Figures 13 and 14, the solution and the relative error are plotted for cases IV and V. 
In both these cases, the agreement with the exact solution is quite good and the relative 
error of the density is approximately .1 percent. In both these cases, the radial variation 
is more pronounced q = 1,2, r^/rq = .65 and k = 1, respectively. Note to more clearly 
demonstrate the radial variation, we have set the amplitude of the vortical waves to be ten 
times what they were in Figures 10-12. 

The amount of transient noise in the solution is very small as is observable from the 
small oscillations in the relative error versus iteration plots. The effect of these transients 
on the computed solution is not discernible in the physical solution. 
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3.5.2 Acoustic Disturbances in a Swirling Mean Flow 

We test the nonreflecting boundary conditions for a variety of swirling flows. The propa- 
gation of an acoustic disturbance is computed and compared with its eigenvalue solution. 
As in the uniform flow case, no reflected disturbances should occur and we examine this by 
computing the transmission coefficients at the downstream boundary. In the ideal case, the 
transmission coefficient of the incident acoustic wave is one at the downstream boundary 
and the reflection coefficients are zero. The inflow specification is given for several differ- 
ent mean flows for the parameters shown in table 2. At the inflow boundary, an incident 


Parameters 

Case I 

Case II 

Case III 

Case IV 

Case V 

WT m / Cm 

24.18 

24.18 

24.18 

24.18 

24.18 

m 

-12 

-12 

-12 

-12 

-12 

Mr 

.5 

.5 

.5 

.5 

.5 

rh/r t 

.65 

.65 

.65 

.65 

.65 

Mxm 

.5 

.5 

.5 

.5 

.5 


.4 

.4 

.4 

.4 

.4 


Table 2: Selected parameters for computing the propagation of sound in a swirling mean 
flow. 

acoustic disturbance with pressure field of the form, 

P = P 1 = Pin(r)e i(me+ktx0) (64) 

is imposed. In Figures (16-18), we present the converged numerical solution for the prop- 
agation of an acoustic disturbance in a mean flow with free vortex swirl. The grid used 
is 80 x 25 x 21 in the x — 9 — r directions, respectively. Figure (16) shows the relative 
error of the calculation versus time iteration and indicates that the converged numerical 
solution lies within half a percent of the exact solution. In Figures (17-18) the numerical 
solution is indicated by the solid line, while the exact solution is denoted by the dotted 
line. Figure (17) shows the solution as a function of axial distance at the intersection of 
the mid-passage and mid-radius planes. Figure (18) shows the solution variation from hub 
to tip. Note that the incident disturbance corresponds to the first radial mode since the 
radial variation of the solution does not cross the axis corresponding to zero pressure from 
hub to tip. The agreement between the numerical and exact solution is good as the solid 
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lines and dotted lines lie on top of one another. Moreover, the agreement remains good 
even as we approach the downstream boundary. In order to obtain a quantitative measure, 
we compute the transmission coefficients at the outflow boundary for the amplitudes of the 
two largest modes, a_i 2 ,i and a_i 2 , 2 * We find the incident mode, ft_ 12 ,i, has an amplitude 
ft— 12,1 = .999. This indicates very little reflection at the downstream boundary. We also 
find that ft_i2,2 = .003, suggesting that the small error in the numerical solution results in 
weak scattering of the incident mode into another mode. 

In Figures (19-21), we present the converged numerical solution for the propagation of 
an acoustic disturbance in a mean flow with rigid body swirl. In this case, the mean flow 
is not potential as in the case of free vortex swirl. Figure (19) shows the relative error of 
the calculation versus time iteration. The numerical solution converges to within half a 
percent of the exact solution given by the normal mode analysis. The numerical solution is 
denoted by the solid line while the exact solution denoted by the dotted line. Figure (20) 
shows the axial variation of the acoustic wave and Figure (21) gives the radial variation 
of the acoustic wave. The agreement shows that the numerical scheme converges to the 
correct solution. The outflow boundary transmission coefficient boundary is ft_i 2 ,i = .999 
indicating very minor reflection at the downstream boundary. Again weak scattering into 
the second radial mode occurs at the boundary. 

Figure (22) shows the axial variation of the converged numerical solution for a case with 
mean swirl that is a combination of free vortex and rigid body swirl. The input disturbances 
was the second radial acoustic mode which crosses the zero pressure axis one time. Again 
the agreement between the numerical solution and the exact solution is good. Computing 
the transmission at the outflow boundary, we find ft_i 2,2 = .997. 

The results presented for the acoustic propagation show that the boundary conditions 
and the numerical interior domain solver converge to the physical solution. These results 
are further verified by the transmission coefficients for the various mean flow distributions 
considered. In the next section, we test the boundary conditions on vortical excitations in 
swirling flow. We also examine the effect that different inflow modeling has on the wake 
evolution and hence on the aerodynamic response for wake blade/row interaction. 


NASA/CR— 2003-2123 1 1 


25 



3.5.3 Vortical Disturbances in a Nonuniform Mean Flow 


In this section, we prescibe vortical disturbances at the inlet and examine the sensitivity of 
the wake evolution to the prescribed inlet conditions. A mean flow distribution represen- 
tative of an engine design is shown in Figure (23). The solution indicates that downstream 
of the fan, the mean flow contains radial variations in both swirl and axial velocity. The 
hub-tip ratio for the case is .58. A disturbance was imposed at the inlet with a reduced 
frequency corresponding to two times blade passing frequency, u — 43. 


We consider inlet excitations where u x ,ug are the same, u r = 0 and we examine the 
effect that different pressure distributions have on the evolution of the vortically-dominated 
wave. The radial pressure distribution takes the form, 


p = p m sin[tn 


(r - r h ) 
(r t - r h ) 


(65) 


where p m = (+.5, 0, -.5). The solid line corresponds to p m positive, the dashed line cor- 
responds to p m zero and the dot-dashed line corresponds to p m negative. Note that the 
case p m = 0, u r = 0 represents the ’pure convection’ approximation which does not satify 
the Euler equations for nonuniform flows. In Figures (24-27), we plot the axial evolution 
of the magnitudes of the pressure and velocity field respectively. We see that the pressure 
magnitude drops very quickly to a small value and does not differ greatly from the zero 
pressure case. The axial and swirl velocities show differences on the order of 20-30 percent 
in their magnitude as they evolve downstream. The large drop in pressure near the exit to 
a small value is consistent with the high frequency analysis which showed that the pressure 
field associated with the vortically-dominated disturbances is small relative to their velocity 
field. These results suggest that the evolution of the upwash component of the velocity is 
only moderately sensitive to the inlet pressure distribution. 


The radial velocity perturbation is not negligible upstream of an exit guide vane. How- 
ever, this information is often neglected from experimetal measurements of the flow field 
between the fan and the stator. To assess the sensitivity of the upwash evolution with 
different inlet condition, we consider different radial velocity distributions with u x ,ug the 
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same, and p = 0. The radial velocity distribution takes the form, 


Ur — U m r 


sin[7ri 


(r - rj . ) 

in - r h ) 


( 66 ) 


where u mr = (+.5,0, -.5). The solid line corresponds to u mr positive, the dashed line u m g 
zero and the dot-dashed line u mr negative. In Figures (27-28), we plot the evolution of the 
magnitudes of the pressure and velocity field, respectively, as they propagate downstream. 
Figure (28), shows that the pressure magnitude drops very quickly to a small value for all 
three cases. This is expected given the high frequency and circumferential mode order of 
the case. However, unlike the previous cases, the magnitude of the pressure is markedly 
different depending on the radial velocity distribution. The radial, swirl and axial velocity 
components are order one and they show differences in magnitude of nearly one hundred 
percent one-half of a mean radius downstream of the inlet. This shows that the solution 
evolution is very sensitive to the prescribed inlet radial velocity distribution. 


3.6 Conclusion 


Inflow and outflow conditions were formulated and implemented for nonuniform mean flows. 
It was observed that the numerical solutions were sensitive to the boundary conditions and 
the details of the inflow excitation. For problems with acoustic waves, the boundary con- 
ditions require accurate solutions for the normal modes in a duct. For vortical excitations, 
the inflow disturbance must be accurately and completely represented. All three compo- 
nents of the wake harmonics must be determined in order to accurately compute the wake 
upwash component and thus the gust response. 


4 Computing the Intensity of Sound in Nonuniform Flow 

4.1 Introduction 

Often one wishes to determine via experiment or computation the spectral composition 
and intensity of sound propagating in a duct [23]. In order to determine the intensity of 
sound, we must separate from the total unsteady field the acoustically-dominated field and 
determine the energy flux of the acoustic waves. In section 2, we used normal mode analysis 
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to represent the unsteady disturbances in a duct and in section 3 we applied the normal 
mode analysis to separate the vortically-dominated field from the acoustic disturbances. In 
this section, we derive an energy relation for acoustic waves propagating in a nonuniform 
mean flow and we use the normal mode solutions to compute the sound power. 

For small-amplitude disturbances to an isentropic irrotational flow an acoustic energy 
equation can be derived from the Euler equations which depends only on first order quanti- 
ties [24, 25]. The intensity flux which results is conserved. Myers [26], derived the acoustic 
energy relation directly from the general energy equation and showed that the acoustic en- 
ergy represents, to leading order, the energy carried by the unsteady disturbances. In the 
general case of a nonuniform mean flow, the time-averaged intensity flux is not conserved. 

In order to study nonuniform mean flow effects on the propagation of sound energy, 
we start from the acoustic energy equation for small-amplitude disturbances and derive, in 
section 4.2, a simplified relation valid in the high frequency limit relevant to the acoustics 
of fan engines. We note that the nonuniform mean flow, in addition to changing the 
aerodynamic interaction mechanism and the acoustic modes, also affects the the expression 
for the acoustic energy as will be discussed in section 4.2. Moreover, as a result of the 
nonuniformity of the mean flow the duct radial modes are not orthogonal and therefore 
interference between the different radial modes will modify the computed sound intensity. * 
In section 4.3, we compute the acoustic power for a homentropic mean flow with free vortex 
swirl to examine the effect that mean flow swirl has on the propagation of acoustic energy 
in a duct. 


4.2 Sound Power in a Nonuniform Flow 


The transport of energy by small-amplitude disturbances is governed by the conservation 
equation [26] 

dE „ r _ 

~dt + V ‘ 1 ~ ( 67 > 

where the energy density, E, the energy flux vector, /, and source term, G, are defined by, 


E = + (po|“| 2 + p'u • Uo), I = (— + u ■ Uo)(pou + p'U), 

*Po Po 


( 68 ) 
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G = «.[^ 0 x(M , -;V)]+^f^-(/Vp 0 -p'Vpo) + ^-(p , ^-p , ^X-69) 

P 0 zp 0 ot ot 

From equation (12), the mean radial velocity is zero and the term, £-f/ 0 • (p'Vpo - p'^Po) 
vanishes. Note that the expressions for the energy and energy flux are second order expres- 
sions which depend only on first order quantities. As a result, the energy of the unsteady 
disturbances can be computed as a by-product of the linearized Euler equations. Expressing 
the last four terms in G in terms of the entropy, equation (67) becomes 


~ + V ./ = «•[(/„ x (p 0 ? - Cop') + £(p'VS„ - j'Vpo)] + (70) 


1 


,dp\ 


Equation (70) shows that the disturbance field interacts with the mean flow and mean flow 
gradients to transfer energy between the mean flow and the propagating disturbances. 


4.2.1 Acoustic Energy in the High Frequency Limit 

In fan engine aeroacoustic calculations, the frequencies of the noise source are often high. 
For example, in tone noise the frequency is (nBQ) where n is the harmonic index, B is 
the number of blades and ft is the shaft rotation frequency. As a result, the characteristic 
acoustic wavelength is £ = c m /{nBQ.) = r m /(nBM m ) where M m is the rotational Mach 
number at the mean radius of the fan. However, variations in the mean flow are character- 
ized by the mean radius of the duct. In the high frequency limit, C/r m << 1. In this limit, 
the vorticity content of the acoustic waves is small [14]. To examine the relative order of 
the various terms in (70), we introduce the fast variables, 


x = x/e, ( 71 ) 

t = t/c, 


where e = t/r m << 1. Substituting into (9) we obtain, 

|V + e(i?-V)So = 0, (72) 

where = J^ + (7 0 -V and V = eV. Equation (72) shows that to leading order the unsteady 
entropy field is convected with the mean flow, 


Do 

Di 


= 0, 


(73) 
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over the small distance, i. 

Substituting (71) into (70), we obtain to leading order 

“ + (T4) 

Note that £ = c£' is 0(|u|). Equation (74) shows that the acoustic intensity is no longer a 
conserved quantity as a consequence of the unsteady entropy and the interaction between 
the unsteady vorticity and the mean flow. When incident waves are only acoustic, then £' is 
small [14] and from (72) s' is 0(cu|). As a result, the acoustic intensity is locally conserved 
over a distance ,0(C), which is small relative to the mean radius. In the general case, the 
disturbance field contains both acoustic and vortical parts. However, in the high frequency 
limit, the two parts can be decoupled since, as was shown in section 3, the pressure field 
associated with the nearly- con vected modes is small. Thus, taking only acoustic part of 
the disturbance field and substituting into (74) we find that the acoustic energy is locally 
conserved, i. e. 

BE 

[-^ + V • /]4 = o, (75) 

where the subscript A denotes the acoustic part of the unsteady field. 

4.3 The Sound Power in a Homentropic Mean Flow 

In this section, we use the normal mode analysis to define acoustic disturbances in a ho- 
mentropic flow. We then derive the sound power and examine the contributions of various 
terms in the sound power for a case of free vortex swirling flow. When calculating the sound 
energy in an annular duct, we are mainly interested in the mean flux of acoustic energy 
across an axial plane. The axial component of the intensity then reduces to 

I* — (p'/po + M x u x + M s iig)(pou x + M x p'). (76) 

Note the expression for the axial intensity has two terms which contain the swirl velocity. 
Substituting (17) into the momentum equations we obtain expressions for the modal velocity 
field in terms of the pressure and radial velocity, 

„ f ( r \J{k mn x+m0) 

Umrir — u mnjmny ' ) c 


NASA/CR— 2003-2123 1 1 


30 



— flfi 


(*dr ^[xfnyn />o) e i(k mn x+mB) 


(77) 


UmriQ — O'Tnn 


[i\j:( rM s)fmn - (rnp m n)/(rpo)\ j (kmnX + me) 
A-mn 


where 


fmn — 


( 2 «^, + A .4( fc , W ) 




(78) 


Substituting (77) into (76) and integrating across an axial surface we obtain the sound 
power in the duct, 




n=N r n=N r m=Ng 

11 =^{E E E a mn a- mn ,e i{ 

71=1 7l # = l 771 = 1 

H 2 PmnPmn' ^ 2 P mnU mn , fi PO^x u mn x u mn f + P()M s U mnx ^ mn f\dr 1 , 

PqCq r ~ 6 6 




r t/ r m 

r [(l + —£-)PmnK 

rh/rm c 0 


(79) 


where iV r and iV# are the number of radial and circumferential propagating modes, * is 
used to denote the complex conjugate and 3? denotes the real part. Since the mean flow is 
potential and isentropic, the time averaged intensity flux is conserved. 

Note that in a swirling mean flow, the eigenfunctions are not orthogonal. As a result, 
interference occurs between different radial modes. Note that the sound power varies with 
x for the interference modes n ^ n'. Thus, the conservation of sound power requires that 
the sum of the modal powers corresponding to the interference modes must be zero, i. e., 

n=N r n=N r , r ,/r m M 2 M MM 

EE r[(l + ^) Pmn U- mn , + -^p mn p* mnl I ( 80 ) 

n^l rjt l Jr hbm C 0 P0C ° C 0 * 

J rpoM x u mnx u‘ mn i i + poM s u mnx u^ nn iJdr = 0 n ^ n' . 

The purpose of the next subsection is to check this result in a homentropic and potential 
mean flow and to determine which of the terms in (79) are the major contributors to the 
sound power. 


4.3.1 Numerical Results for a Free Vortex Swirling Flow 

In a free vortex swirling flow, acoustic disturbances propagating in a duct are purely po- 
tential and conserve acoustic energy. To examine the major contributions to the sound 
power and study the effect of swirl, we compute the contributions to the modal sound 
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power from each of the terms in equation (79). We consider a free vortex swirling flow in 
an annular duct, as shown in Figure (32). The plot shows that the annulus extends from 
a non-dimensional radius of .8 to 2.0 and the swirl Mach number is inversely proportional 
to the radius. To conserve sound power, the sum of the contribution of each interference 
term inside the integral in equation (80) must be zero so that there is no x-dependence in 
the sound power. Thus, although the eigenfunctions are not orthogonal, the contribution 
of interference terms to the sound power should be negligible. 

In Figure (33), we show the contributions of each of the terms in the integral for the 
downstream propagating modes corresponding to a case with reduced frequency, u = 20 and 
spinning mode order, m = -1. The open circle, + and * correspond to the contributions of 
the first three terms in (80), respectively and the x and diamond correspond to the fourth 
and fifth terms. The abscissa corresponds to the radial mode index n' and the ordinate 
represents the contribution of each term. The different figures correspond to different 
radial mode orders, n for a given spinning mode order, m. The index, n, describes the 
eigenfunction which has n - 1 zero crossings from the hub to tip. In this case, there exist 
6 propagating modes. The contribution of the interference terms are quite small relative 
to the noninterference terms, n = n' . Also the dominant contributions to the power result 
from the contributions of the first, second and fourth terms. For the highest order radial 
mode, n — 6, the first two terms contribute the most to modal power. As the mode order 
decreases, the contribution of the first term increases rapidly from a value of .25 to nearly 
one. The second term does not change appreciably for the various mode orders but the 
contribution of the fourth term becomes significant for the lower order radial modes. For 
this case, with low spinning mode order, m = - 1, the third and fifth terms which contain 
the swirl component of the mean velocity produce negligible contributions to the modal 
power. Summing the noninterference terms of the six modes we find the contribution is of 
the order 1 x 10 -3 as opposed to the noninterference terms which result in 0(1) contributions 
to the power. 

In Figure (34), we consider a case with a reduced frequency of 50 and spinning mode 
order m = -16. In this case, there exist 15 propagating modes which are again ordered 
based on the number of zero crossings that occur between the hub and the tip. We consider 
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the six upstream propagating modes which range from n = 9 to n = 14. Again the first, 
second and fourth terms contribute the most to the sound power. As in the previous 
case, the contribution of the fourth term increases for the lower radial modes. Interference 
between the different radial modes was observed but the magnitude of the interference 
terms were O(10“ 2 ) as opposed to the noninterference contributions which are 0(1). The 
most significant interference contributions occurred from the modes adjacent to the n = n r 
mode. 

Finally, we consider a case with reduced frequency, 24.5 and spinning mode order m = 
—28. In Figure (35), we show the modal power contributions for each term in (79) . In this 
case, the swirl plays a more significant role in the modal power. The third term in (79) 
is most significant for the higher order radial modes. Contributions from the interference 
terms are evident in this case. However, their summed contributions remain small, i. e. 
0 { 10 " 3 ). 

4.4 Conclusions 

An expression for the sound power was derived which is valid for nonuniform flows. For 
high frequency disturbances, the propagation of sound energy can be treated as locally 
conserved. The dominant contributions to the sound power were computed for several un- 
steady acoustic disturbances. These results can be used in conjunction with an aeroacoustic 
calculation to compute the sound power and determine the major modal contributors to 
noise. For the case of a free vortex mean flow swirl, interference effects were negligible. 
However, these results must be examined for more general nonuniform mean flows. 
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Figure 1: The radial profile of the mean velocity and total enthalpy for a plane lying 
between a fan and an exit guide vane. The computational results were obtained using a 
Reynolds averaged Xavier Stokes code and show the strongest radial variations in the tip 
region of the annulus where the tip vortex of the fan lies. 
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Figure 2: Effect of Error in Nonreflecting Boundary Conditions on Numerical solution. Row 
1 shows the global relative error in the numerical solution versus time. Row 2 shows the 
real part of the numerical solution along the one-dimensional domain for increasing errors 
in the boundarv condition, Ak. 
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Figure 3: Plots comparing a theorectical meai 
vortex and rigid body swirl. The mean flow m 
data is represented by the ircles. 
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Figure 4: Comparison of eigenmodes spectra obtained using the spectral method and the 
shooting method, m = —1, u = 16. The plus symbol denotes the initial guess obtained 
by the spectral method and the open circle denotes the final converged solution which is 
obtained by the shooting method. Note that the spectral method shows five upstream and 
downstream propagating modes, while the shooting method indicates that there are only 
four propagating modes. 
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Figure 7: Eigenmode spectra for different mean flow models, m = — 12, uj = 24.2. Note 
that the mean flow modifies the upstream propagating modes more than the downstream 
propagating modes. 
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Figure 8: Eigenfunctions of upstream propagating acoustic modes, m — —12 , cj = 24.2; 
v mean flow data. theoretical mean flow model. Note that the theoretical mean flow 
model provides an accurate representation of the higher order radial modes but does not 
do well for the first two radial modes. 
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First Downstream Acoustic Mode Second Downstream Acoustic Mode Third Downstream Acoustic Mode 



Figure 9: Eigenfunctions of downstream propagating acoustic modes, m = —12, u? — 24.2; 

mean flow data, theoretical mean flow model. The eigenfunction solution using the 
theoretical mean flow model provides a good approximation except for the first radial mode. 
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Figure 10: A superposition of a vortical and an acoustic disturbance is imposed upstream. 
The hub-tip ratio is .65 and the radial phase variation is such that a quarter of a wavelength 
fits between the hub and the tip. The numerical solution of the density and two components 
of the velocity are plotted as solid lines as a function of x , 6 and r at the midplanes of 
the duct respectively. The exact solution is shown by the dashed lines. The global relative 
error of the density is plotted in the bottom right figure as a function of time iteration. 
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I iguie 11. A superposition of a vortical and an acoustic disturbance is imposed upstream. 
I'lie hub-tip ratio is .65 and the radial phase variation is such that a half of a wavelength fits 
between the hub and the tip. The numerical solution of the density and two components 
°f l * u ' velocity are plotted as solid lines as a function of x, 6 and r at the midplanes of 
the duct respectively. The exact solution is shown by the dashed lines. The global relative 
error of the density is plotted in the bottom right figure as a function of time iteration. 
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Figure 12: A superposition of a vortical and an acoustic disturbance is imposed upstream. 
The hub-tip ratio is .65 and the radial phase variation is such that a full wavelength fits 
between the hub and the tip. The numerical solution of the density and two components 
of the velocity are plotted as solid lines as a function of x, 6 and r at the midplanes of 
the duct respectively. The exact solution is shown by the dashed lines. The global relative 
error of the density is plotted in the bottom right figure as a function of time iteration. 
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Figure 13: A superposition of a vortical and an acoustic disturbance is imposed upstream. 
The hub-tip ratio is .65 and the radial phase variation is such that a full wavelength fits 
between the hub and the tip. The amplitude of the acoustic wave is 1.0 The numerical 
solution of the density and two components of the velocity are plotted as solid lines as a 
function of x, 6 and r respectively. The exact solution is shown by the dashed lines. The 
global relative error of the density is plotted in the bottom right figure as a function of time 
iteration. 
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Figure 14: A superposition of a vortical and an acoustic disturbance is imposed upstream 
The hub-tip ratio is .65 and the radial phase variation is such that two full wavelengths fits 
between the hub and the tip. The numerical solution of the density and two components 
of the velocity are plotted as solid lines as a function of x, 6 and r at the midplanes of 
the duct respectively. The exact solution is shown by the dashed lines. The global relative 
error of the density is plotted in the bottom right figure as a function of time iteration. 
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Figure 15: The velocity field of the two times blade passing frequency harmonic result- 
ing from a wake shed from an upstream fan. The computational results were obtained 
using a Reynolds averaged Navier-Stokes code. Note that the velocity field is fully three- 
dimensional and that the axial, circumferential and radial velocity components are all the 
same order of magnitude. 
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Figure 16: The convergence history for the propagation of an acoustic wave in a free vortex 
swirling flow. The reduced frequency, spinning mode order, and hub-tip ratio for the case 
are 24.18, -12. .65. respectively. The convergence parameter, E, is the relative error in the 
Filler calculation where the ’exact’ solution is the eigensolution computed from the normal 
mode analysis. The free vortex mean flow is chosen such that the swirl Mach number at the 
mean radius is .4. 


NAS A/CR— 2003-2 12311 


52 







Figure 17: The axial variation of the converged Euler solution is compared with the eigen- 
solution obtained from the normal mode analysis. The numerical solution is denoted by 
the solid line and the ’exact’ solution is denoted by the dashed line. The good agreement 
between the numerical and analytical solution even at the exit of the domain show the 
accuracy of the numerical scheme and the nonreflecting boundary conditions. The Euler 
solution contains eighty points in the axial direction. 
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Figure 18: The radial variation of the Euler solution compared with the eigensolution 
obtained from the normal mode analysis. The numerical solution is denoted by the solid 
line and the ’exact’ solution is denoted by the dashed line. The Euler solution contains 
twenty points from hub to tip. 
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Figure 19: The convergence history for the propagation of an acoustic wave in a rigid body 
swirling flow. The convergence parameter, E , is the relative error in the Euler calculation 
where the ‘exact’ solution is the eigensolution computed from the normal mode analysis 
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Figure 20: The axial variation of the Euler solution compared with the eigensolution ob- 
tained from the normal mode analysis. The Euler solution is denoted by the solid line and 
the ’exact 1 eigensolution is denoted by the dashed line. The good agreement between the 
numerical and analytical solution even at the exit of the domain show the accuracy of the 
numerical scheme and the nonreflecting boundary conditions. The Euler solution contains 
eighty points in the axial direction. 
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Figure 21: The radial variation of the Euler solution compared with the eigensolution 
obtained from the normal mode analysis. The numerical solution is denoted by the solid 
line and the ’exact’ eigensolution is denoted by the dashed line. The Euler solution contains 
twenty points from hub to tip. 
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Figure 22: The axial variation of the Euler solution compared with the eigensolution ob- 
tained from the normal mode analysis. The Euler solution is denoted by the solid line and 
the exact eigensolution is denoted by the dashed line. The good agreement between the 
numerical and analytical solution even at the exit of the domain show the accuracy of the 
numerical scheme and the nonreflecting boundary conditions. The Euler solution contains 
eighty points in the axial direction. 
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Meanflow 



Figure 23: The swirl and axial mean flow distribution representative of the flow between 
fan and a fan-exit guide vane for a jet engine. 
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Pressure Evolution 



Fi gure 24: The evolution of the pressure amplitude for three different gust inflow conditions 
described by u x = / x (r), u e = / 2 O), u r = 0 and p = p m sin[7ri ^~ r r ^ ] where /i(r),/ 2 (r) are 
given and p m = +.5,0. -.5. 
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?o: The evolution of the radial velocity amplitude for three diffe 
ns described by u x = fi(r),ug = f 2 (r),u r = 0 and p = /> m sin[; 
(r) are given and p m = +.5,0, -.5. 
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Figure 27: The evolution of the swirl velocity amplitude for three different gust inflow 
conditions described by u x = fi(r),u e = / 2 (r), tz r = 0 and p = p m sin[ 7 ri^^^] where 
are given and p m = +.5,0, -.5. 
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Pressure Evolution 



Figure 28: The evolution of the pressure amplitude for three different gust inflow conditions 
described by u x = fi(r),ug = /2(r),u r == u mr sin[7rf ^73^ ] and p = 0 where /i(r),/ 2 (r) 
are given and u mr = +.5,0, -.5. 
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Figure 29: The evolution of the radial velocity amplitude for three different gust infl 
conditions described by u x = fi(r),ug = / 2 (c), u r == u mr sin[rf ] anc j p = q w hi 

/ 1 ( r). f 2 ( r) are given and u mr — +.5, 0, -.5. 
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Figure 32: Free vortex swirl distribution for a duct with hub-tip ratio of .4 
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Figure 33: Contributions of each term in (79) to the sound power of the downstream 
propagating duct modes. The reduced frequency is lj - 20 and m = -1. The abscissa is 
the radial mode index, n' and each figure corresponds to the contribution of each radial 
mode n. The top row going from left to right correspond to n = 6, 5, 4 in descending order 
and the bottom row going from left to right n = 3,2, 1. 
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Figure 31: Contributions of each term in (79) to the sound power of the downstream 
propagating duct modes. The reduced frequency is u = 50 and m = -16. The abscissa 
is the radial mode index, n' and each figure corresponds to the contribution of each radial 
mode n. The top row going from left to right correspond to n = 16, 15, 14 in descending 
order and the bottom row going from left to right n — 13, 12, 11. 
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Figure 35: Contributions of each term in. (79) to the sound power of the downstream 
propagating duct modes. The reduced frequency is u = 24.5 and m = -28. The abscissa 
is the radial mode index, n' and each figure corresponds to the contribution of each radial 
mode n. The top row going from left to right correspond to n = 6,5,4 in descending order 
and the bottom row going from left to right n = 3,2, 1. 
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